A multi-species model with interconversion, chipping and injection 
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Motivated by the phenomenology of transport through the Golgi apparatus of ceUs, we study 
a multi-species model with boundary injection of one species of particle, interconversion between 
the different species of particle, and driven diffusive movement of particles through the system by 
chipping of a single particle from a site. The model is analysed in one dimension using equations for 
particle currents. It is found that, depending on the rates of various processes and the asymmetry in 
the hopping, the system may exist either in a steady phase, in which the average mass at each site 
attains a time-independent value, or in a "growing" phase, in which the total mass grows indefinitely 
in time, even in a finite system. The growing phases have interesting spatial structure. In particular, 
we find phases in which some spatial regions of the system have a constant average mass, while other 
regions show unbounded growth. 
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I. INTRODUCTION 

Living cells possess multiple trafficking path- 
ways, which have in common, a regulated ffux of 
cargo molecules, such as proteins and lipids, mov- 
ing through and processed within, organized com- 
partments or organelles. For instance, the secre- 
tory pathway consists of molecules that move from 
the endoplasmic reticulum (ER) to the cell surface 
via an organelle system called the Golgi appara- 
tus |T]. The Golgi apparatus itself consists of dis- 
tinct sub-compartments known as cisternae. Pro- 
teins and lipids arrive from the ER, enter through 
one face of the Golgi and undergo several chemical 
reactions (processing); the modified products then 
leave through the other face to mainly target the 
cell surface. 

There has been much discussion about what 
drives the flux of molecules through the Golgi [2]. 
The vesicular transport model envisages that the 
cisternae are stable structures with fixed enzymatic 
composition. Molecules shuttle from one cisterna to 
the next in small sacs called vesicles, and get chem- 
ically modified by the resident enzymes. The cis- 
ternal maturation model, on the other hand, consid- 
ers the cisternae to be transient structures that are 
formed by fusion of incoming vesicles. In this model, 
it is the cisternae that progress through the Golgi ap- 
paratus, carrying the biomolecules with them. Spe- 
cific enzymes get attached to a cisterna in different 
stages of its progression, and modify its contents. 
The final cisterna eventually breaks up, releasing 
processed biomolecules. Independent studies on a 



variety of cells provide evidence for both these pos- 
sibilities [3H5]. 

Not only the cargo molecules, but indeed the 
molecules that form the Golgi organelle themselves, 
must also be trafficked along the same route and by 
the same driving forces. This invites the following 
question: How does the Golgi organelle form in the 
first place, i.e., how does one obtain stable structures 
(cisternae) given the rules of molecular trafficking, 
which broadly, may be described as: (i) localised 
injection of 'particles', i.e. of the vesicles contain- 
ing unprocessed biomolecules (ii) transformation of 
particles from one species to the other, i.e. chemi- 
cal processing of the biomolecules by enzymes (iii) 
transport of particles either by chipping (breaking 
off) of a single particle (corresponding to vesicle 
movement) or through movement of bigger aggre- 
gates (corresponding to cisternal progression). The 
aim of this paper is to construct a statistical model 
incorporating these elementary processes, and use 
this to quantitatively address questions of structure 
formation and the nature of the states at long times. 

With this motivation, we define the following 
multi-species model. Particles of species A are in- 
jected into a one-dimensional (ID) lattice at one 
boundary. Particles of different species B, C... (or 
more generally, of all types A,B,C etc.) leave from 
the other boundary. This happens by allowing A 
particles to convert to B particles (and vice versa), 
B particles to C and so on. There is no restriction on 
the number of particles of any species a site can hold. 
The hopping of particles from one site to another can 
either occur collectively via movement of the whole 
stack or one at a time by chipping of a single particle. 



Chipping refers to a single particle breaking off from 
a stack and hopping to a neighbouring site. The 
hopping probability may be the same to the left and 
right (diffusive) or different (driven diffusive, due to 
an existing chemical or electrical gradient). When 
a particle (or a collection of particles) hops on to a 
stack, it merges with the particles already resident 
on that site. Thus, stacks constantly gain and lose 
particles. 

This is a generalization of a well studied model 
of aggregation and chipping [5', T . Earlier studies 
dealt with a closed system, with a single species of 
particle. The present generalization deals with an 
open system, with injection of particles at one end, 
and interconversion from one species of particle to 
another. Interestingly, we find that new sorts of 
phases can arise in some limits. 

The parameter space is large. Thus, it is useful 
to begin with the study of the 'chipping only' model, 
where there is no movement of stacks as a whole. In 
the remainder of this paper, we will study the model 
in this limit, with chipping and interconversion rates 
taken to be constants, independent of the number of 
A or B or C... particles on the site. With this as- 
sumption of constant rates, we find that for some 
rates, the system fails to achieve steady state in the 
sense that unbounded growth of mass occurs. Inter- 
estingly, even in these growing states, the particle 
currents at each site are stationary (time indepen- 
dent) after sufhciently long times. The indefinitely 
growing average mass at a given site arises simply 
because the particle currents, though stationary, are 
not balanced at that site [5]. Thus, we call such a 
state quasi-stationary. 

Although we have defined the model for an 
arbitrary number of species, from now on we will 
focus primarily on the two-species model. The 
multi-species model is a simple generalization of the 
two-species case and shows qualitatively similar be- 
haviour, as discussed in Sec. |VI| 

The rest of the paper is organised as follows. 
Section |ll] defines the model precisely, highlights 
some connections with other models, and briefly dis- 
cusses the main results of the paper. In Sec. |III[ we 
analyse the behaviour of the first site in detail. In 
Sec. |IV) we study the case of fully asymmetric hop- 
ping for a ID lattice and solve the current equations 
at each site. In Sec. [V] a continuum limit is taken 
and results obtained for general asymmetry. In Sec. 
|VI[ we briefly discuss the generalization of the model 
to three species. Section |VII| contains a discussion 
of some issues and extensions. 



II. THE MODEL AND RESULTS 
A. Model 

The model is defined on a ID lattice. At any time 
i, a lattice site i has mf particles of type A and mf 
particles of type B. We start with an empty lattice 
with L sites at t = 0. At each instant, a site i is 
chosen at random and one of the following moves 
(illustrated in Fig. [ij occurs in an infinitesimal time 
interval At: 

i. Injection: If the site picked is site 1, then an A 
particle is injected into it from the left bound- 
ary with probability a At: ra^ -t- 1. 

ii. Interconversion: With probability uAt, one of 
the A particles residing on the selected site 
converts to type B: if mf > 0, then mf — > 
mf — 1 and mf m,f + 1. With proba- 
bility vAt one of the B particles residing on 
the site converts to type A: if mf > then 
mf — >■ mf — 1 and mf — > mf + 1. 

iii. Chipping: With probability p_R At {p^At), one 
of the A particles chips off the site and hops to 
the right (left) neighbouring site, thereby de- 
creasing the mass of A on the site by 1 and 
increasing the mass of A on the neighbour- 
ing site by 1; if mf > then mf — > mf — 1 
and mf^i — >■ mf^^ + 1. With probabihty ^At 
(g^At), one of the B particles chips off the site 
and hops to the right (left): if mf > then 
mf ^ mf - 1, mf±i ^ mf^^ + 1. 

iv. Stack Movement: With probability D^At (or 
D^At), the entire mass (i.e. all A and B par- 
ticles together) present at the site hops to the 
right (left) neighbouring site: mf — > 0, mf — 
0, mf?ti ^ mfj^i + mf, mf^^ mf^^ + mf. 

At the last site an A (B) particle can hop out of the 
system at a rate PR{qR), i.e. with the same hopping 
rate as in the bulk. Once particles hop rightwards 
from the last site L, they cannot return to the sys- 
tem. Thus, site L + 1 acts as a sink. 

In general, the drive i.e. the asymmetry in 
the rates of the three hopping processes - hop- 
ping of stacks, hopping of A particles after chip- 
ping and hopping of B particles after chipping, can 
be different. However, for simplicity, we will take 
the asymmetry to be the same for all three pro- 
cesses, i.e. Pr/pl = Qr/ql = Dr/Dl. Then we 
can parametrize the rates as follows: pr, — jp, 
Pi = (1 - l)p; IR = 79, 9L = (1 - l)q; Dr = 71}, 
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FIG. 1. Illustration of the allowed moves with A (white) and B (black) particles: Injection of A particles at rate a at 
site 1. A ^ B and B ^ A conversions at rates u and v respectively (see site 4). An A particle can hop rightwards 
(e.g. from site 5 to 6) at rate pr and leftwards (to site 4) at rate pl- Similarly, a B particle (from site 7), can hop 
to sites 8 or 6 at rates qn and respectively. Stack movement corresponds to all particles on a given site (here site 
2) collectively hopping to the left or right neighbour (site 1 or 3) at rate Dl or Dr respectively. 



Dl = (1 —j)D. Here, 7 is the measure of the asym- 
metry in hopping and takes on values between and 
1. 

The moves (i) to (iv) above define a general 
model in which particles can move via both stack 
movement and chipping. However, as indicated in 
section |Tj in this paper we will analyse the model 
in the absence of stack movement i.e. in the limit 
D = 0. 

Several features of our model, including bound- 
ary effects, mass movement and coalescence, and 
two-species coupling have been studied separately in 
a number of different contexts. Macroscopic aggre- 
gates are known to form in models such as the the 
zero range process (ZRP) ^9J and the aggregation- 
chipping model O [7] . In these models, the closed 
system in steady state can enter a condensate phase 
in which a finite fraction of the total mass resides on 
a single site. A quasi-stationary state characterised 
by growth of total mass in the system occurs in ag- 
gregation models such as the Takayasu model 
in which particles are injected at every site. The 
in-out model [HJ [11], which extends the Takayasu 
model to allow for evaporation of unit masses, ex- 
hibits a phase transition from a growing state to a 
steady state. Such a phase transition is also seen in 
our model. However, the phases in our model are 
not spatially uniform as injection occurs only at the 
boundary. 

Open boundaries with injection can induce in- 
teresting features such as boundary-driven phase 
transitions in the asymmetric simple exclusion pro- 
cess (ASEP) [iT. In the ZRP with open boundaries 
|13j . steady state, in the sense of constant average 
mass, is not reached for strong boundary injection; 
indefinitely growing aggregates may form on one 



or both boundaries, while bulk sites attain steady 
state or show slower growth, depending on the in- 
jection, hopping rates, and asymmetry. Our model, 
in which injection takes place at the boundary, also 
shows growing phases; somewhat counter intuitively, 
a steady state with constant mass is reached close 
to the boundary, while there is unbounded growth 
of the mass in the bulk. 



Coupling between different species of particle 
may or may not involve interconversion. For exam- 
ple, in the two-species ZRP [T31[TS], condensates of 
one or both species may form above critical densi- 
ties but the number of particles of each species in the 
system does not change. New effects are seen when 
interconversion (or more generally, non- conservation 
of particle number) occurs. Multi-species models 
that allow for inter-conversion have been studied ear- 
lier in the context of transport in single P^^E(7| and 
multiple pTH21] channels |25|. These models differ 
from the present work in that they restrict the occu- 
pancy to at most one particle per site in each chan- 
nel. Generically in an open boundary channel, if 
the number of particles is not conserved in the bulk 
(e.g. if switching between channels or deposition and 
evaporation occurs), then for low switching rates, 
the system may develop a localised shock separating 
low density and high density regions in one of the 
channels [231 121] • Our model too shows co-existence 
of two phases of one species (a steady phase and a 
growing phase) in different spatial regions separated 
by a domain wall with a tunable location. 
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FIG. 2. Schematic depiction of a phase with pile-ups of B for (a) fully asymmetric (7 — 1) (b) fully symmetric 
(7 = 1/2) hopping. Dashed lines show (m^) as a function of x at two time instants ti (small dashes) and t2 (big 
dashes) respectively such that ^ ti < t2- The two lines overlap in regions of steady state (i.e. (m^) remains 
constant in time); the line lies above the ti line in regions of B pile-ups {{mB) increases with time). Insets; 
Occupation probability as a function of x. While = 1 in regions of B pile-ups, < 1 in regions of steady 
state. For 7 = 1/2, the region of pile-ups (xi < a; < X2) exists between regions of steady state near each boundary; 
for 7 = 1, the entire region to the right {x > Xi) of the steady state region has B pile-ups. 



B. Results 

Our main results may be summarised as follows: 
The model is parametrized by the rates a,p,q,u,v 
and the asymmetry parameter 7, where 7 = 1/2 
corresponds to symmetric and 7 = 1 to fully asym- 
metric rightward hopping. In different regions of pa- 
rameter space, different kinds of behaviour are seen: 

i. For some values of parameters, the system ex- 
ists in a steady phase. Occupation probabili- 
ties of A and B at a site i, i.e. the probability 
that the site has at least one A (B) particle are 
defined as 
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where Pi (m^ 



,t)dt is the probability of 



finding A particles and B particles on 
site i between time t and t + dt. In the steady 
phase, and sf at each site, reach a time 
independent value less than 1. The total mass 
in the system (and on each site) also reaches a 
finite time-independent value. 



occupancy s 
(see Sec. llTl) 



For other values of parameters, the system ex- 
ists in a quasi-stationary state characterised 
by stationary (time-independent) currents but 
non- stationary i.e. indefinitely growing total 
mass. This happens even in a finite system. 
We will refer to unbounded growth of mass as 
formation of pile-ups. More precisely, if the av- 
erage A (B) mass at a site grows indefinitely, 
we refer to this as an A (B) pile-up at the site. 
At the site with an A (B) pile-up, the mean 
"^{s^) approaches 1 at long times 
This may be treated as a func- 
tional definition of a pile-up. A pile- up is to 
be distinguished from a condensate of the sort 
found in the Zero Range Process (ZRP) [3]: A 
condensate at a site contains a finite fraction of 
the total number of particles in the system but 
the condensate mass does not grow in time. 

iii. In the growing phase, a system may have pile- 
ups of A, of B or of both species. For the situa- 
tion with injection of only A particles, as stud- 
ied in this paper, an A pile-up can be found 
only on site 1 (except in the special case 7 = 1, 
g = 0), i.e. it is a boundary phenomenon. B 
pile-ups, on the other hand, are found in the 
bulk. In a phase with B pile-ups, an interest- 
ing situation can occur where the lattice has a 
region in steady state ( with finite mean mass ) 
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co-existing with a region of B pile-ups (where 
the mean mass grows indefinitely). For fully 
asymmetric hopping, the region close to the 
left boundary is in steady state while pile-ups 
occur in the entire region to the right (shown 
schematically in Fig 2(a) ). For the symmetric 



or partially asymmetric case, there are two re- 
gions in steady state close to either boundary, 
with a region of B pile-ups between these (Fig. 



2(b)) 



5^=a/w_ 



4=^ 



In the rest of the paper, we will adopt the fol- 
lowing terminology. If the average mass of both A 
and B reaches a finite time-independent value, we 
will refer to the system as being in a steady phase. 
If one region of the system attains constant average 
mass and another region shows unbounded growth, 
we will refer to the two regions as a region in steady 
state and a region with pile-ups respectively. The 
system as a whole will be referred to as being in a 
growing phase. 




(a) 




(b) 



FIG. 3. Schematic diagram showing variation of the par- 
ticle currents , and j with distance x from left 
boundary for 7 = 1 for a) a system in steady phase b) a 
system with a steady state region and a B pile-up region. 
In a region of steady state, total particle current j does 
not change with x. In a region of pile-ups, j{x) decreases 
with X 



FIG. 4. In the absence of interconversion (it = n = 0), 
the occupation probability profiles and Sq (dashed 
lines) are flat, being a/w and respectively. When u 7^ 
0, V 7^ 0, there is a net interconversion from A to B 
at each site, pushing down from the value and 
pushing up from the value Sq . Here and (solid 
lines) are occupation probabilities with interconversion. 



Steady and growing phases can also be distin- 
guished by the spatial profile of the currents. Let 
jf, jf and ji ~ if + denote respectively the 
net A particle current, the B particle current and 
the total particle current in the ith bond of the lat- 
tice i.e. from site i to i + 1. Figs. 3(a) and |3(b)] 
show schematically how, for 7 = 1, these three cur- 
rents vary across the lattice in a steady phase and 
a growing phase (with B pile-ups) respectively. In 
the steady phase, although and vary with the 
distance from the origin a;, they do so in such a way 
that the total particle current j remains constant 
across the lattice. In contrast, in the growing phase, 
the total particle current j also varies with x in the 
pile- up region. In fact, ji^i must be greater than 
ji for a pile- up to exist at site i, implying that j 
decreases with i (or x) in a pile-up region. 

In order to explain the spatial structure of the 
various phases, it is useful to understand what causes 
the spatial variation in j"^ and or equivalently, s"^ 
and . For simplicity, let us consider purely right- 
ward hopping (7 = 1) in a system of size L. In the 
absence of interconversion, the particle current 
is the same in each bond, and so is . Thus, s^ 



and s are constant across the lattice, being a/w 
and respectively (dotted lines in Fig. |4]). Intercon 
version changes and 
s^ (solid lines in Fig. 



" and consequently s 
in the following way. 



^ and 

_ The 

particle current into site 1 (jo) is entirely of type A. 
At site 1, some of the A particles convert to type B, 
resulting in a small probability that site 1 also has B 
particles. Hence, ji is different in composition, hav- 
ing a small B component as well. At site 2, there is 
still a much larger concentration of A particles than 
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B particles, and hence a net conversion of A to B 
particles. Thus, the composition of j2 shifts further 
towards B. This process continues until we reach a 
site at which there is no net conversion. Beyond this 
site and do not change with x. If intercon- 
versions are faster than hopping (i.e. u,v > p,q), 
then at a given site, the relative amounts of A and 
B can alter significantly through interconversion be- 
fore hopping occurs. In this limit, and change 
sharply over a few sites and reach their asymptotic 
values very close to the left boundary. In the other 
limit {u,v <^ P,q), and vary slowly over a 
large part of the lattice. In particular, if u and v 
are not 0(1) but 0(1/L), then the spatial variation 
in and becomes a bulk rather than boundary 
effect. 

The right boundary is an absorbing boundary. 
The degree of asymmetry in the hopping determines 
how far into the lattice the effect of the absorbing 
boundary condition propagates. For fully asymmet- 
ric dynamics (7 = 1), the right boundary has no ef- 
fect on any other site while for 7 slightly less than 1, 
the effect extends only into a very small region of the 
lattice close to the boundary. For 7 = 1/2, the in- 
formation about the right boundary propagates into 
the entire lattice. Thus, the asymmetry in hopping 
rates and the rate of interconversion relative to the 
rate of hopping together determine how and 
or s"^ and vary across the lattice. 



III. ANALYSIS OF THE FIRST SITE 

For fully asymmetric (rightward) hopping (7 — 
1), a site is not affected by any sites to the right 
of it. Thus, the first site can be analysed indepen- 
dently of the rest of the lattice; this analysis yields 
considerable insight into the full problem. We find 
four possible behaviours: steady state; pile-up of A 
but not B; pile-up of B but not A; pile-ups of both 
A and B. 

The master equation for the probabil ity d istribu- 
tion P{m^ ,m^) (as defined in section II B I at the 
first site can be written as follows: 



For m,-^ > 1 and > 1 
9P(m'4,m^^ 



dt 



+ qP{m^, + 1) + aP{m^ - 1, m^) 
+ uP{m^ + 1, - 1) + vP{m^ - 1, 
-{a + u + v+p + q)P{m^, m^) 



1) 

(2a) 



^=pP(l,m^)+gP(0,m^ + l) ^2b) 

+ uP{l, -l)-{a + v + q)P{0, m^) 



dP{m^,0) 
Ft 



= pP(m^ + 1,0)+ gP(TO^, 1) 



+ aP(TO^ -1,0) + wP(m^ - 1, 1) 
- {a + u + p)P{m'^,0) 



(2c) 



5P(0,0) 
dt 



pP(l,0) +gP(0,l) -aP(0,0) (2d) 



To understand the behaviour of the system, we 
do not solve the master equation but work with the 
time-evolution equations for (m"^) and (m^): 



d{m^) 
dt 



ps 



dt 



us 



qs 



(3a) 



(3b) 



Here (m"^) and (m^) denote the average masses of 
A and B respectively on site 1 averaged over time 
histories and s^ and s^ the occupation probabili- 
ties of A and B (defined in Eq. ([T])); all of these 
are functions of time. Evidently, Eq. ([s]) is just a 
continuity equation for the average masses. 

The time-evolution of the first site can also be 
mapped to a random walk in two dimensions. The 
two-dimensional (2D) space here is the m'^ — 
space. The random walk starts at the origin at t = 0, 
moves in both positive and negative directions along 
the axis (with rates a and p respectively) but 
moves only in the negative direction along the 
axis (with a rate q) as there is no injection of B par- 
ticles. Diagonal moves, corresponding to intercon- 
version, are also allowed but only in the direction of 
increasing and decreasing or vice versa. The 

and axes are reflecting boundaries. A possi- 
ble random walk, with all allowed moves in m^ — m^ 
space, is illustrated in Fig. |5(a)| Clearly, the ran- 
dom walk cannot be decomposed into independent 
walks along the m"* and axes because of the 
diagonal moves. 

The random walk picture is useful in the fol- 
lowing way. To begin with, let us consider a simple 
ID random walk in the positive half of the a;-axis 
with a reflecting barrier at the origin. This walk has 
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(a) 




(b) 



FIG. 5. Time evolution on the first site can be mapped to random walk in 2D — space, (a) A possible 
random walk in — space. Allowed moves are enumerated in the upper right corner of the figure along with 
the corresponding rates at which they occur. 

(b) Trails in — space corresponding to various scenarios. Trail in solid line shows steady state, trail (in dashed 
line) directed along axes shows pile-up of A, trail (in dashed line) along axes shows pile-up of B, fourth trail 
(dot-dash) shows pile-up of both A and B. 
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FIG. 6. V ~ u phase diagram for the first site when A and B hopping rates are equal (a) p = q,p/a = 0.2 (b) 
p — q,p/a — 0.8 



a drift velocity c and a diffusion constant D. Af- 
ter a length ^ £'/|c|, drift becomes more important 
than diffusion. If c > 0, then beyond this length, 
the mean displacement will grow linearly with time, 
while if c < 0, the mean displacement reaches a con- 



stant value ^ D/\c\. If c = 0, then the motion is 
always diffusive, and the mean displacement grows 
as ^/i. The 2-d walk is more complicated because 
motion along and axes is coupled due to the 
diagonal moves (corresponding to interconversion) . 
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The diagonal moves effectively make c and D in the 
direction dependent on the average displacement 
along the direction. If the effective c in the 
direction is greater than or equal to zero, then we 
get a pile-up of A on the first site (with c = im- 
plying growth of the pile-up and c > 0, growth 
that is linear in time). Also, as expected from the 
theory of random walks, the fluctuation about the 
mean mass, i.e. {(m^Y) — {m^)'^ grows as ^/t (for 
c > 0). If c < 0, then the mean mass of A and the 
mean deviation in mass of A reach time indepen- 
dent values. In terms of the occupation probability, 
if < 1, then must be finite and if s"* = 1, 
is infinite at long times (i.e. pile-ups occur). To see 
this, note that the return from to'^ = M to = 
requires M steps, implying that P{m^ = 0) falls as 
g-oAf^ which is non- vanishing only for finite M. So, 
if = 1 — P{m^ = 0) is less than 1 at long times, 
< > must be finite; if is equal to 1 at long 
times, < > is infinite. 

From the above discussion it is clear that de- 
pending on the rates, the effective c in either direc- 
tion (i.e. along and axes) can be negative 
or non-negative, resulting in four possibilities alto- 
gether. Figure |5(b)| shows schematically, possible 
random walks corresponding to the four scenarios - 
steady state, pile-up of A, pile-up of B and pile-up 
of both A and B. 

In terms of the rates, the conditions for each of 
the four scenarios can be derived from Eq. ([3|. For 
example, if there is a pile-up of B but not of A on 
site 1, then at long times, d{m^)/dt > with — 1 
and d{m-^)Jdt — with < 1. Substituting these 
into Eq. (p| gives the two inequalities in parameter 
space that must be satisfied in this case. Appendix 
[A| lists the necessary conditions in parameter space 
for each of the four cases. 

These conditions can be visuahsed better by 
considering a reduced parameter space. We look at 
a case where both species hop with the same rate, 
i.e. p = q. Then the total mass evolves according 
to: 



d((m^ 



dt 



(4) 



Clearly, if a < p, then s , s < 1 and the system will 
always attain steady state. If a > 2p, the system can 
never attain steady state and must show pile-ups (as 
d{{m^) + {m^)) / dt > even when and take on 
the maximum possible value i.e. 1). In this case, if 
u i.e. the interconversion rates are similar, both 
A and B pile up. However, if one of the conversion 
rates is much higher than the other, then the species 
which converts fast is able to stabilize, so that only 



one of the species piles up. Figure 6(a) shows the 
V — u phase diagram for the case a > 2p. li p < a < 
2p, then and can adjust themselves such that 
steady state is attained. However, as before, if one 
of the species converts much faster than the other, 
it drives the la tter o ut of steady state, causing it to 
pile up. Figure [6 (b) | shows the v — u phase diagram 
for the case p < a < 2p. 



IV. ASYMMETRIC HOPPING ON A ID 
LATTICE 

In this section, we analyse the fully asymmetric 
case (7 = 1). In this limit, the current equations in- 
volving and have a simple solution. Also, this 
case illustrates in a clear way, how interconversions 
drive spatial gradients in A and B particle currents. 

Current equations similar to Eq. ([S]) can be 
written for each site in the system: 



d < mf > 
dt 



US: 



ps^ 



d < mf > 
Jt 



and for site 1, as before 

d < > 



dt 



BAA 

a + vsi — usi — psi 



d < 



dt 



US-i — f Si — 



qsf 



(5a) 

(5b) 
(6a) 
(6b) 



The above equations are exact, being just the con- 
tinuity equations for the average masses. They can 
also be obtained by writing a master equation for the 
total probability V{{mf}, {nif}), multiplying it by 
and averaging over all configurations. 



A. Steady Phase: 



For the steady phase {d{m^)/dt — d{m^)/dt = 
at all sites), equations ([5| and ^ give recursion 
relations expressing sf and sf in terms of sf_i and 
•sf-i- Iterating these recursion relations gives the 
following steady phase occupation probabilities: 



uq + vp 



uq + vp 



[1 - A'] 



(7) 
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FIG. 7. Occupation probabilities sf and sf as a function 
of site number i in a steady phase for 7 = 1, a = 1, 
u = 0.5, V = 0.25, p = 2, q — 1. Data points are 
obtained from Monte Carlo simulations, and solid lines 
are plots of the analytical results in Eq. ([7|. and 
approach their asymptotic values over a length Id 
defined in Eq. ([8|. 



where A — {pq)/{uq + pv + pq). 

Figure [t] shows the spatial profile of sf and sf 
in a steady phase, as obtained from Monte Carlo 
simulations, along with the analytical results, for a 
certain choice of parameters. 

From Eq. Q, the following can be deduced: 

i. In a steady phase, occupation probabilities sf, 
sf must be less than 1 for all i. This is ensured 
if the maxima of the two expressions in ([7]) are 
less than 1, i.e. a{v + q)/{uq + pv + pq) < 1 
(for s^) and au/{uq+pv) < 1 (for s^). Thus, 
having rates which satisfy these two inequali- 
ties simultaneously is a necessary and sufficient 
condition for the system to be in steady phase. 

ii. Away from the left boundary, sf and sf ap- 
proach asymptotic values that satisfy the rela- 
tion sf / sf = v/u which corresponds to zero 
interconversion current at site i. 

iii. The characteristic length [Id] associated with 
this (exponential) approach is given by 



In 



(8) 



which is {v/q + u/p)~^ to first order in inter- 
conversion rates. This first order term has the 
following interpretation. 1/p is the time over 
which a single hopping event of an A particle 
takes place at each site in the lattice, u is the 
rate of conversion for an A particle. Thus, u/p 
is the probability that a conversion event (of 



an A particle) takes place in this time. The 
term v/q has a similar interpretation. Thus, 
{v/q + u/p) is an estimate of the average num- 
ber of interconversions taking place in the sys- 
tem in a time interval over which or before 
hopping occurs. As reflected in the expression 
for I o (and discussed in section II B I , when 
this number is large, s^ and s^ approach their 
asymptotic values over just a few sites from the 
left boundary. Conversely, for interconversion 
rates that are proportional to 1/L, the length 
Id is of order L and the gradients in s"*, ex- 
tend over the bulk of the lattice. (0(1/L) in- 
terconversion rates are also studied in [231124].) 



B. Growing Phases: 

To study the structure of the growing phases 
in the system, we use Fig. |8] which is a diagram- 
matic version of Eqs. ([5| and ([6]). Each site can 
be thought of as having A and B compartments. 
The figure shows particle currents between the vari- 
ous compartments, with vertical arrows representing 
intra-site interconversion currents and horizontal ar- 
rows the inter-site hopping currents. The currents 
into and out of any compartment can either balance 
each other (resulting in a steady state for the com- 
partment) or the in current can be more than the 
out current (resulting in a pile- up). Evidently, the 
out current cannot exceed the in current at any site 
as this would eventually lead to negative mass. This 
simple observation can be used to establish the fol- 
lowing result: 

In both steady and growing phases, if q ^ 0, then 

\fi. Consequently, 



< si 



and sf_^-^ > s. 



us 



i+1 ' 



usf^i < usf — vsf , i.e. interconversion cur- 



rent keeps decreasing as we move right, asymptoti- 
cally approaching zero. 

To prove this result, we note that ii q 7^ 0, then 
at site 1, the interconversion current, being the only 
'in' current for compartment B, is always positive, 
i.e. from A to B. This implies that sf < sf and 
sf > sf . We can prove this by ruling out the other 
three possibilities which are: 

(i) sf>sf,sf>sf: This would imply that either 
the net current out of site 2 is more than the net 
current into it which is not possible; or if sf = sf 
and sf = sf , then usf — vsf = usf — vsf > 0, 
i.e. there is a net interconversion current out of A. 
Since the in and out hopping currents for A are equal 
{sf — sf), this means there is a net out current from 
compartment A which is not possible. 
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FIG. 8. Each site can be thought of as having A and 
B compartments. Exchange of particles takes place be- 
tween compartments at a give site, i.e. through intercon- 
version currents (represented by vertical arrows). Ex- 
change of particles takes place between sites through 
hopping currents (represented by horizontal arrows). 



(ii) S2 < s^, S2 < sf : Since and sf can at 
most be equal to 1, this implies that s^,s^ < 1, 
i.e. site 2 is in steady state and currents for both 
compartments balance. To satisfy current balance 
for A with < s^, interconversion current must be 
from A to B (see Fig. [s]) but for current balance for 
B with S2 < sf , it must be in the opposite direction. 
This leads to a contradiction. 

(iii) S2 > s^, $2 < sf: This implies that ms^ — 
vsf > us^ — vsf > 0. However (from Fig. [sj, to 
have sf > sf, the interconversion current has to be 
from B to A or at least has to be zero. This leads 
to an inconsistency. 

Thus, the only possibility is that sf < sf and sf > 
sf . In the same way, the proof can be extended 
to all subsequent sites. It may be noted that we 
have not assumed the existence of steady state i.e. 
d{m^'^) /dt ~ 0, anywhere in this argument. Thus, 
the inequalities hold for the growing phase also. The 
(7 = case has many sub cases and can be analysed 
in a similar way. In the rest of the discussion, it is 
assumed that q 0- 

The above inequalities for the occupation prob- 
abilities have the following implications: 



i. Even if a pile-up of A exists on site 1 {sf = 1), 
no subsequent site can have a pile-up of A (for 
q =^ 0) a.a sf < sf Wi > 1. Thus, an A pile- up 
would be found only at site 1, if at all. 

ii. The occupation probability s^ increases from 
left to right. If s^ touches 1 at some site k, 
then sf = 1 Vz > fc as sf > sf V« > k. Thus, 
there are B pile-ups at all sites right of k. 

It also follows from Eq. ([5 1 that when B pile- 
ups exist, then at long times, the mass of B grows 
linearly with time at all sites with pile-ups (except at 
site fc, the earliest site with a B pile-up, where it may 
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FIG. 9. Results of Monte Carlo simulations for a = 1.3, 
u — 1, V — 0.33, p — 0.86, 5=1. Snapshots of the 
lattice showing profile at ti = 200000 (squares) and 
t2 = 300000 (circles) {t in units of 500 MC Steps). Inset: 
Semi-log plot of mf (t2) — mf (ti) vs i. Solid line shows 
the analytical prediction. 



grow as \/t). Moreover, if ™f_(_]^ ^ ct, then nif^2 ^ 
{p/{p + u)) ct, mf^3 ^ (p / {p + u))'^ ct etc. Thus, 
after sufficiently long times, a snapshot of the lattice 
at any instant would show an exponential decay in 
space in the mass profile of B, site A: -I- 1 onwards, 
with a decay length given by [In ({p + u) /p)] . The 
location of the site k depends on the parameters and 
can even be tuned to be at site 1. Left of site fc, the 
system exists in steady state (except for a possible A 
pile-up on site 1) while right of it, the system exists 
in a growing phase. 

The system was also studied numerically by do- 
ing Monte Carlo simulations for small lattices, typ- 
ically L — 20. Since later sites cannot affect earlier 
sites, this is the same as doing a simulation with 
large L but keeping track of only the first few sites. 
Numerical results are presented for a choice of pa- 
rameters for which the theory predicts the existence 
of steady state at the first four sites, and pile-ups of 
B site 5 onwards. Figure |9] shows snapshots of the 
B profile of the lattice at two times ti — 200, 000 
and t2 — 300, 000 (where t is in units of 500 Monte 
Carlo (MC) steps). In agreement with the analyti- 
cal prediction, upto site 4, reaches a finite time- 
independent value, while growth of starts oc- 
curring from site 5, is fastest at site 6 and becomes 
slower as we go further from the boundary. The in- 
set of Fig. [9] shows m^{t2) — m^(ti) as a function of 
site number in a semi-log plot. The solid line is the 
analytically expected exponential decay. There is 
good agreement between the solid line and the data 
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points until site 10 but not after that. For sites after 
site 10, the long time behaviour, i.e. linear growth 
does not set in during the run time of the simulation 
(1.5 X 10^ MC steps). Such large relaxation times 
were the primary reason for choosing small system 
sizes in simulations. 

Thus, to sum up, for only rightward hopping 
(7 = 1), with injection of A particles at site 1, de- 
creases monotonically while increases monotoni- 
cally with site number. The system may attain sta- 
tionarity or it may exist in a growing phase charac- 
terised by 'pile-ups'. An A pile-up can occur only at 
the first (boundary) site (for q ^ 0). B pile-ups occur 
in the bulk; in fact, it is possible to have steady state 
existing near the right boundary and pile-ups of B 
further on (as in Fig[9|. By making interconversion 
much slower than hopping {u/w + v/q ^ 0{l/L)), 
both the steady state and the pile-up region can be 
made macroscopically large {0{L)). If a site has a 
pile-up of B, then all sites right of it also have B 
pile-ups. These pile-ups grow &s ^ t. In the pile- 
up region, at any instant, the spatial profile of 
shows an exponential decay. 

When q = 0, i.e. B particles cannot move, 
there are two important differences. Steady state 
on site 1 implies a steady state behaviour for the 
whole system, and it is not possible to have a region 
of steady state followed by a region of pile-ups. Sec- 
ondly, g = is the only surface in parameter space 
on which pile-ups of A need not be localised at site 1. 
If there exists a pile-up of A on site 1, then as long 
as M < w, there are pile-ups of A on all subsequent 
sites. These pile-ups grow as ^/t with time. 



values between and 1. In the continuum limit, 

1 P)aA-B 1 f)2„A,B 



L dx 



2X2 



A. 7 = 1 

First we would like to check whether the steady 
phase profile given by Eq. Q is recovered in the 
continuum limit. Substituting from Eq. ([9| into 
Eqs. ([5]) and (|6]), setting time derivatives equal to 
zero and retaining only leading order terms in l/L, 
we get: 



ds^ 
dx 



vs^ 



ds^ ~ A 
q—— + us^ 
dx 



Hs^ ^ 



vs^ = 



(10a) 



(10b) 



along with the boundary conditions s^{x = 0) = a/p 
and s^{x = 0) = 0. Here, u and v refer to the 
rescaled rates u — uL and v — vL. In Eq. ( 10 ) 



the full derivative with respect to x has been used 
instead of the partial derivative as the steady state 
occupation probabilities are functions of x only. The 



solution of Eq. ( 10 1 is given by 



\x) 



av 



uq + vp 



uq 
— I 
vp 



■ exp 



\x) 



uq + vp 



1 — exp 



(11a) 



(lib) 



CONTINUUM LIMIT 



For fully asymmetric hopping (7 = 1), steady 
phase occupation probabilities were obtained by 
solving Eqs. ([s]), ([6]). For any other 7 however, a 
given site is affected by both its left and right neigh- 
bours and the resulting time evolution equations for 
< mf > and < mf >, (similar to Eqs. (|6|) 
yield complicated recursion relations for the occupa- 
tion probabilities which are difficult to solve. Thus, 
to analyse the system for any general asymmetry, we 
assume the total number of sites L in the system to 
be large, and take a continuum limit for the lattice, 
thereby going from the integer valued site number 
i to the real valued position co-ordinate x, defined 
as X — i/L. By this definition, x can take on real 



with / = (i/q + u/p) = {l/L){v/q + u/p) . A 
comparison with the decay length computed for the 
discrete lattice (in Eq. ([s])) shows that apart from 
the length rescaling factor 1/L, I is just the first 
order term in the expansion of the exact decay length 
Id in terms of {v/q + u/p). This first order term 
arises because in retaining only the leading order 
term in Eq. ^ , we are making the assumption that 
s^'^{x) are slowly varying functions of x which in 
turn requires that interconversions are much slower 
than hopping, i.e. {v/q + u/p) « 1. 



B. 



7 = 1/2 



The other limiting case is 7 = 1/2 or perfectly 
symmetric hopping. As before, current equations 
similar to Eqs. ^ and ([6| can be written for the 
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discrete lattice, which on taking the continuum Umit, 
yield the following two coupled differential equa- 
tions: 



d < m'^{x,t) > 
dt 



2 



2„B 



d < m"{x,t) > _ q d\s 



dt 



2 dx^ 



(12a) 



(12b) 



Here the rescaled interconversion rates are u = uL^ 
and V = vL'^ . The boundary conditions are given by 
5-^(0) = 2a/p, s-^(O) = and 5-^(1) = s^{l) = 0. 
The boundary condition at x = 1 just reflects the 
fact that the right boundary acts as a sink, as par- 
ticles cannot return to the system once they leave 
from site L. Equation ( 12 ) can be solved for steady 



phase by setting time derivatives equal to zero, tak- 
ing second derivative of Eq. (12a) with respect to 
X and eliminating the terms by using Eq. (12b) 
and the original Eq. (12a). This gives an equation 



which is the second derivative with respect to a; of a 
second order differential equation in s"^. The solu- 
tion of this equation gives the following steady phase 
spatial profiles of and : 
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10(a) 



(1 _ a;) + [ sinh (v^(l - x)) } 
vp y smh{y/r]) j 
(13a) 



(1 



/sinh {y/rj{l - x)) 



sinh (\/^) 



u/p) 



(13b) 



Figure 10(a) shows the spatial profile of s and 



s in steady phase for a certain choice of parameters. 



Clearly, steady phase profiles (Eq. (13)) obtained 
from the continuum approximation tally well with 

the Monte Carlo results. 

The spatial profiles of and in Fig. 10(. 
give an indication of the nature of the growing 
phases of the system. As in the asymmetric case 
with q ^ 0, since s"^ is always a decreasing function 
of x, it follows that an A pile-up can exist only on 
the first site. If an A pile- up exists at site 1, then 
we use the boundary condition s'^(a; = 1/L) = 1 in- 
stead of the boundary condition for {x — 0) in Eq. 
(12) to solve for the steady state of the rest of the 
10(a) we can also deduce that if 
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FIG. 10. Occupation probabilities af and sf as a func- 
tion of site number i 

(a) in a steady phase (with a = 1.3, u = 0.02, v — 0.01, 
p = 2.2, g = 0.5). 

(b) in a growing phase (with a — 1.0, u — 0.0136, v — 0, 
p = 2.2, q = 0.3). 

Data points obtained from Monte Carlo simulations seem 
to agree well with analytical predictions (solid lines). 



the maximum in the s^{x) spatial profile crosses 1, 
the system cannot exist in a steady phase but must 
show pile-ups of B. Moreover, these B pile-ups oc- 
cur in some central region of the lattice flanked by 
regions of steady state on either side. To test this 
prediction, simulations were done with parameters 
for which s^{x), as obtained from Eq. (13), crosses 
1 in some range of x values. Figure 10(b) shows 



system. From Fig. 



the occupation probability profiles for such a choice 
of parameters. As expected, the region of pile-ups 
(with = 1) exists in the middle of two steady 
state regions. Appendix [B] presents the details of 
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how we analytically calculate the spatial boundaries 
of the pile-up region as well as the occupation prob- 
ability profiles which are shown as solid lines in Fig. 
|10(b)[ The main assumption involved in this calcu- 
lation is that ds^ /dx vanishes at the boundaries of 
the pile-up region. 



C. 1/2<7<1 

Away from the limits of pure drift (7 = 1) or 
pure diffusion (7 = 1/2), current equations similar 
to Eq. ( 12 ) can be written with partially asymmetric 
hopping rates (1/2 < 7 < 1). These equations also 
have terms with first derivatives of and with 



respect to x, and can be solved in the same man- 
ner as Eq. (12), giving the following expressions for 



the steady phase occupation probabilities s (x) and 
s^{x): 



s^ix) 



j{uq + vp) 



1 - exp(-2|(l - x)) 



1 - exp(-20 



+ — exp(C) 
vp 



sinh \y\[f]~+^{l - xyj ^ 
sinh i\fi)^^j J 



(14a) 



7(ug -I- vp) 



■ exp(0 



1 - exp(-2|(l - x)) 
1 - exp(-2|) 



sinh y + C 



(14b) 



where u 
f =(27- 



2{v/q + u/p) and 



l)L 



Equation (14 1 is also valid for < 7 < 1/2. 
However, if particles have a net drift back towards 
the source, only a small region close to a; = is 
occupied. Thus, we confine our analysis to the more 
interesting case 1/2 < 7 < 1. Strictly at the two 
limits, 7 = 1/2 and 7 = 1, Eq. (14) is not valid. 
However, taking the limit 7 — > 1/2, or equivalently 



FIG. 

P = 
fj = 
with 



1 

0.9 
0.8 
0.7 
0.6 
0.5 
0.4 
0.3 
0.2 

11. 




8 10 12 14 16 18 20 
site number i 



st and sf for a = 1.2, u = 0.2, v = 0.1, 
2.2, q = 0.7 and 7 — 0.7. For these paramet ers , 
103 and ^ = 8.4. Solid lines are plots of Eq. ([l4| 
X replaced by i/L. 



(13). At 7 = 1, a second order description cannot 



be used as there is only one boundary condition for 



(or s 



The spatial variation of s (x) and s (x) is 
governed by the rescaled variables fj and ^. If 
fj ~ 0[L^) and ^ ^ ^{L), this variation is con- 
fined to the boundary regions, eventually becoming 
a discontinuity in the limit L— >oo. If77^0(l) 
and ^ ~ C'(l), the gradients due to the two bound- 
aries extend into the whole lattice and together, de- 
termine the behaviour of the system in a complex 
way. For large but not infinite values of fj and ^ 
{fl ~ O(IOO), i - 0(10)), gradients near the two 
boundaries can be treated independently and yet are 
not mere bou nda ry discontinuities. This situation is 
shown in Fig. 
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For such values of fj and ^ and suf- 
ficiently large L, two distinct length scales emerge, 
one characterizing the variation near x = and the 
other near x — 1. These length scales are, respec- 
tively: 



2^ 



(15) 



The length li is the distance from a; = beyond 
which A and B are effectively decoupled (no net in- 
terconversion) . li is small for fast interconversions 
(high values o f 77), which has been rationalized ear- 
lier in section 



II B 



li also decreases as ^ decreases, 
i.e. with decrease in hopping asymmetry. 

The effect of the absorbing right boundary (at 
a; = 1) extends over a length which is given by 12- 
This can be seen in the following way. Since A 
and B are effectively decoupled beyond h, the spa- 



for Eq. (14) gives back the profiles of Eq. tial variation of s (x) (and similarly s (x)) beyond 
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this length is described by a drift diffusion equation 
{p/2{d^s^/dx^) - piids'^/dx) ~ 0) along with the 
boundary condition s^{x = 1) = 0. The diffusion 
length associated with this equation is [2^]^^ which 
is just 12- 

As h and I2 increase, the variations near the two 
boundaries can no longer be treated independently 
and the middle re gion of nearly constant (x) and 



s (x) seen in Fig. 11 disappears. Perfectly symmet- 
ric hopping, for which I2 diverges, is an extreme case 
of this. 
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VI. MULTI-SPECIES MODEL 

The two-species model can be extended quite 
easily to more species. For example, consider the 
three-species generalization. Particles of type A are 
injected at the left boundary. An A particle can 
convert to B and vice versa; in addition, B particles 
can convert to C (at rate k) and C to B (at rate I). 
Note that direct interconversion between A and C 
is not allowed. Also, as before, a single particle of 
any species can chip off a site (with rates p, q, r for 
species A,B,C respectively) and move, in general, in 
a driven diffusive way. 

The model can be analysed in exactly the same 
way as the two-species model by writing equations 
involving particle currents. The analysis of the 
three-species model yields the following results: 

i. While decreases monotonically as a func- 
tion of a;, and s'-^ both show qualitatively 
the same behaviour of the spatial profile, i.e. 
both increase monotonically with x in the 
driven case and have a peak in the diffusive 
case. This indicates that in the three-species 
model, depending on the rates, pile-ups of B 
or C or of both can occur in the bulk of the 
lattice and a pile-up of A only on the first site. 

ii. Interestingly, although interconversion takes 
place in the order A ^ B ^ C ^ pile-ups need 
not appear in the same order in space. For ex- 
ample, C pile-ups may be found to the left of 
B pile-ups (i.e. closer to the left boundary) as 
shown in Fig [12] 

Thus, in the multi-species model with sequen- 
tial interconversion, in different regions of parameter 
space, the system may have pile-ups of one, a few or 
all (except A) species of particle. Pile-ups of dif- 
ferent species may occur in different regions of the 
system. 



FIG. 12. Results of Monte Carlo simulations for a phase 
with pile-ups of both B and C (7 = 1). Parameters: 
a = 5, ^i; = 6, (7 = 2, r- = 0.3, u = 0.7, v = 0.2, k = 0.65, 
I — 0.62. Data points show mf and mp at two time 
intervals ti = 5 x 10^ and t2 — IQ^ ■ Pile-ups of C occur 
earlier in the lattice (site 13 onwards) than pile-ups of B 
(site 15 onwards). Inset: sf, sf , at long times. 



VII. CONCLUSION 

In this paper, motivated by the phenomenol- 
ogy of traffic in the Golgi apparatus of the cell, we 
have studied a stochastic two-species model with 
boundary injection of type A, interconversion be- 
tween types A and B, and transport of both species 
through the lattice (in general, in a driven diffusive 
way) by chipping of one particle at a time. We found 
that depending on the rates of various processes, the 
system may either eventually attain steady state or 
show unbounded growth of mass at all times. Pile- 
ups (as defined in the text) may be composed of A 
particles or B particles or both. Unlike in transla- 
tionally invariant systems, the phases in our model 
have interesting spatial structure. Generically, in 
growing phases, a part of the lattice attains steady 
state while other regions show unbounded growth of 
mass (pile-ups). 

Below, we comment on some possible extensions 
of the model: 

i. If stack movement also occurs (as in the more 
general model defined in section II A I, there 



are no growing phases and the average mass 
in the system always attains a constant finite 
value. This may be argued as follows. If jin 
and jout are the in and out currents at a site, 
then the average mass m at the site grows as 
dm/dt = jin —jout- In the 'chipping only' case, 
iout = ps'^ (assuming only one type of parti- 
cle). If jin > p, then in and out currents can- 
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not balance, and m grows indefinitely. On the 
other hand, with stack movement, jout which 
is given hy ps^ -\-Dm^ ^ can always balance 
resulting in zero growth rate or constant mass. 
However, if the rate of stack movement is mass 
dependent i.e. D{m) cx m~", then jout which 
is now ps^ + D{fn-^Y~°' will be bounded for 
a > 1. In this case, pile-ups can occur and the 
behaviour of the system would be expected to 
be similar to that of the 'chipping only' model. 

ii. One extension of this model to higher dimen- 
sions is trivial. If injection takes place at a sur- 
face perpendicular to one of the spatial direc- 
tions, and we assume periodic boundary con- 
ditions in other directions, the model is still 
effectively one dimensional as there are no net 
currents between sites in the transverse di- 
rection. Generalizations to higher dimensions 
with other boundary conditions can be more 
complex. 

iii. This model can also be extended quite easily 
to include injection of B particles at the left 
boundary, at a rate h. Now, depending on the 
value of the ratio auq/bpv, the system either 
has net interconversion from A to B (and be- 
haviour similar to what has been discussed in 
the paper) or net interconversion from B to A 
(resulting in A pile-ups in the bulk etc.) or 
no net interconversion, i.e. no effective cou- 
pling between the two species. At the right 
boundary, introducing exit rates that are dif- 
ferent from the chipping rates can also change 
the behaviour of the system. However, this 
change is only a boundary effect unless the exit 
rate is 0{1/L) times the chipping rate. If the 
exit rate is small in the sense defined above, 
then s"^ and vary with x differently from 



Eqs. (13). However, growing phases with B 



pile-ups in the bulk are still found. 

Finally, we comment on the possible relevance 
to the biological system. The model system we have 
studied is a simple one, but it shares the following 
qualitative features with the Golgi apparatus. First, 
there is a gradation of sizes of stacks from the source 
to the sink. In particular, small stacks are found 
close to the source while large, though growing, ag- 
gregates occur in regions farther away. Secondly, 
there is a gradation of 'chemical' species across the 
system, with different spatial regions having aggre- 
gates that may be predominantly of type B and/or 
type C etc. These aggregates, however, show un- 
bounded growth, which is not realistic in the context 



of the Golgi, which consists of large but bounded and 
discrete structures. As discussed above, unbounded 
growth is eliminated once we allow for movement of 
whole stacks, as in the more general model defined 
in Sec III Al and discussed in this section. 



In the present model, we have assumed that 
mass transfer occurs either by chipping of a single 
particle or by movement of the whole stack. Other 
transport processes could involve the fragmentation 
and movement of chunks of intermediate sizes. The 
study of these different mass transport processes in 
systems with boundary injection is an interesting di- 
rection for future study. 
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Appendix A: Possible behaviours on the first 
site: 



As indicated in Section pTH conditions for each of 
the four scenarios for the first site can be derived 
from Eq. ([s]). They are as follows: 



i. Steady state: If 



a{v+q) 



< 1 and 



< 



qp-\-qu-\-pv qp-\-qu-\-pv 

1, then steady state is reached with steady 
state probabilities: 



qp-{-qu-\-pv 



qp-i-qu-i-pv 

ii. Pile-up of B: If > 1 and < 1, 

qp+qu-\-pv — u-\-p ' 

then (m^) keeps growing while (m^) reaches 



a constant value. 
At long times: 



dt 



ua—vp 
u+p 



gA ^ a±v 
u-\-p ' 



iii. Pile-up of A: If > 1 and ^ < 1, 

qp+qu+pv — v+q ' 

then {m^) keeps growing while (m^) reaches 
a constant value. 



= 1 and 



At long times: 



v+q ' 



— 1 and 



dt 



a — p 



uq 
v+q 



iv. Pile-ups of both A and B: If ^ > 1 and 

^ u+p — 

> 1, then both (m"^) and (m^) keep 
growing. At long times: 
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and 



s = s 

d{r. 



~ dt 



V — U — p 



dt 



Appendix B: Method for calculating s''^{x) and 
s^{x) for 7 = 1/2 when pile-ups of B occur: 

Let xi and X2 divide the lattice into three regions 
such that in: 

Region I (0 < a; < xi): Steady state {s^{x) < 1. 
s^ix) < 1); 

Region II {xi < x < X2)- Steady state of A, pile-ups 
of B {s^{x) < 1, s^{x) = 1); 

Region III {x2 < x < 1): Steady state {s^{x) < 1, 
s^ix) < 1). 

In regions I and III, since steady state exists, 
s^{x) and s^{x) are obtained by setting the L.H.S. 
equal to in both of Eqs. ( |12a ) and ( 12b ) and solv- 
ing them with the following boundary conditions for: 
Region I: s^(0) = 2a/p, s^(0) = 0, s^(xi) = 1, 
ds^/dx\ =0 



Region III: s^{x2) = 1, ds^/dx\ 
s^(l) = 



= 0, 



(1) 



In region II, substituting s^{x) = 1 into Eq. 
(12a), and setting the left-hand-side equal to zero, 



we get s (x) upto two unknowns. We can find 



s (xi) and s {X2) from the solutions in regions I 
Since 



and III respectively. Since s (x) must be a con- 
tinuous function, these provide the boundary con- 
ditions for the solution of s^{x) in region II. Thus, 
we have s^{x) and s^(x) for all three regions upto 
the two unknowns xi and X2. To determine xi and 
X2, we make the assumption that ds^/dx is con- 
tinuous at xi and X2, i.e. ds^/dx\ = ds"^/dx\ 



and ds^ ldx\^^ = ds^ Idx T ' ' which yields two tran- 

'1X2 '1X2 ■' 

scendental equations in xi and X2. For a given set 
of parameters, these can be solved numerically to 
get xi and X2, Since we have expressions for s^{x) 
and s^(x) in terms of xi and X2, this also gives the 
occupation probabilities in all three regions. 
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This calculation involves two ad hoc assump- 
tions, namely, that ds^ /dx vanishes and ds^ /dx is 
continuous at the two boundaries of the pile-up re- 
gion. These assumptions give results which match 
well with numerical data (see Fig. |10(b)[ ), but they 
need to be investigated and adequately justified. 
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